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Abstract 

Effects of feedbacks on self-organization phenomena in networks 
of diffusively coupled bistable elements are investigated. For regular 
trees, an approximate analytical theory for localized stationary pat- 
terns under application of global feedbacks is constructed. Using it, 
properties of such patterns in different parts of the parameter space 
are discussed. Numerical investigations are performed for large ran- 
dom Erdos-Renyi and scale-free networks. In both kinds of systems, 
localized stationary activation patterns have been observed. The ac- 
tive nodes in such a pattern form a subnetwork, whose size decreases 
as the feedback intensity is increased. For strong feedbacks, active 
subnetworks are organized as trees. Additionally, local feedbacks af- 
fecting only the nodes with high degrees (i.e. hubs) or the periphery 
nodes are considered. 

1 Introduction 

Nonequilibrium pattern formation has been extensively studied in distributed 
active media [l}j3]. In excitable or bistable media, stationary localized pat- 
terns can develop. They represent domains with high activator density sur- 
rounded by large regions where the density of the activator is much lower. 
Such structures have been theoretically considered EJIH] and experimentally 
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observed in chemical reactions [6}[7] and in semiconductors |8l|9] . They have 



also been studied under the global feedback control 10 11 



Control of nonequilibrium patterns, as well as their purposeful design, is 
an essential issue that has attracted much attention 12 -14 . Global or local 



feedback control schemes serve as a standard method used for this purpose. 
Typically, global feedback requires a common control signal, generated by 
the entire system and applied back to all its elements. In contrast, local 
feedback control does not directly affect the entire system, but is applied 
only to some of its elements. Various feedback schemes have been used 
in theoretical studies fTolfTolfTol and in the experiments [l~7T[20l, either for 



|10|[15|[16] and in the experiments 
stabilizing existing unstable patterns or for inducing new kinds of patterns 
that do not existing in the absence of feedback. 

Self-organization phenomena, such as epidemic spreading 21 25 , clus- 
tering 26 and synchronization 



or traveling and pinned fronts 30 



27 



28,29 



of oscillators, Turing patterns 
have been studied in reaction-diffusion 
systems organized in complex networks. Moreover, some effects of control by 
global feedbacks have been previously investigated for networks systems. It 
has been demonstrated, for example, that turbulence in oscillator networks 
can be suppressed [3l] and hysteresis of Turing network patterns can be 
prevented [32] when such feedbacks are applied. 

As we show in this Letter, feedback control may also suppress spreading 
of activation in networks of bistable elements, and leads to the emergence 
of localized stationary patterns which resemble stationary spots in continu- 
ous media. The approximate analytical theory of such phenomena could be 
constructed for regular trees and systematic numerical simulations for ran- 
dom Erdos-Renyi (ER) and scale-free (SF) networks were undertaken. Both 
global and local feedbacks, applied to a subset of network nodes, were con- 
sidered. Properties of developing stationary patterns, depending on control 
conditions, have been explored. 



2 Bistable systems on networks 

Classical continuous one-component bistable media are described by u(x, t) = 
f(u, h)+DV 2 u('x, t) , where w(x, t) is the activator density, D is the activation 
diffusion constant and function f(u, h) specifies the local bistable dynamics. 
For example, this function can be chosen as f(u, h) = u(h — u)(u — 1) , where 
the parameter h determines the activation threshold and plays an important 
role in front propagation. 

In network-organized systems, the activator species occupies the nodes 
of a network and can be diffusively transported over network links to other 
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nodes. The architecture of a network is described in terms of its adjacency 
matrix T, whose elements are = 1, if there is a link connecting the nodes 
i and j = 1, N), and 7y = otherwise. We will consider processes in 
non-directed networks, where the adjacency matrix is symmetric (7y = Tji). 
Generally, the network analog of one-component bistable system is given by 

N 

iii = f(ui, hi) + D ^(TijUj - TjiUi) , (1) 
i=i 

where Ui is the amount of activator in network node % and f(ui, hi) describes 
the local bistable dynamics of the activator. The last term in eq. takes 
into account diffusive coupling between the nodes and the coefficient D char- 
acterizes the rate of diffusive transport of the activator over the network links. 
An essential property of a node i is its degree fcj (the number of connections) 
given by h = Y.j T ji- 

We introduce the global feedback through the parameter h, i. e. as 

h = h + f i,(S-So), (2) 

where n > is the intensity of the feedback, S = Yl!j=i u j * s the total ac- 
tivation of the network and So is a parameter defining the size of localized 
patterns. In our simulations, it was taken equal to the number of the nodes 
which were initially activated. Hence, the threshold h depends now on the to- 
tal activation. It increases when more nodes are activated, so that a negative 
feedback is realized. 



3 Regular tree networks 

We first investigate the system (jl]) without feedback control, where spreading 
or retreating activation fronts can be found. Fronts can also get pinned, 
thus forming stationary patterns. Previously, we have derived the pinning 



condition for regular trees depending on the parameters k and D 30 , but 
now we would also need to know how the pinning behavior depends on the 
threshold h. We analyze this behavior and determine how it is affected by 
the feedback. 

Let us consider the system on a regular tree with the branching factor 
k — 1. In such a tree, all nodes, lying at the same distance I from the origin, 
can be grouped into a single shell. Suppose that we have taken a node which 
belongs to the shell I. This node should be diffusively coupled to k — 1 nodes 
in the next shell / + 1 and to just one node in the previous shell / — 1, as 
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Figure 1: (Color online) Stationary patterns corresponding to pinned fronts 
of type I (a) and type II (b) in a regular tree with the branching factor 
k — 1 = 2. By grouping the nodes with the same distance from the root into 
a single shell, the two types of fronts can be displayed as shown on the right 
side. The parameters are (a) h = 0.35 and (b) h = 0.65; D = 0.02 



described by the equation 

in = f(u h hi) + D(m-i - ui) + D(k - l)Om - ui) , (3) 

where ui is the density of the activator in the shell I. Note that for k = 2, 
eq. Q3J) corresponds to a one-dimensional chain of coupled bistable elements, 
where pinned fronts have been previously investigated j33|[34). However, the 
respective approximate analytical theory for the trees (k > 2) has been only 



recently developed 30 



Generally, both the activation front propagating from the root to the 
periphery of a tree (type I) and in the opposite direction, i.e. towards the tree 
root (type IT), are possible. First, we derive approximate pinning conditions 
for these two types of fronts, depending on the branching factor k — 1 and 



the parameter h. If diffusion is weak enough (cf. 30,34]), a pinned front can 
be found by setting iii = in eq. ([3]), so that we get f(ui,hi) + D(ui~\ — 
Ui) + D(k — l)(ui + \ — Ui) =0. Suppose that a front of type I is pinned 
and located at the shell / = m. The front is so sharp that the nodes in 
the lower I < m shells are all approximately in the active state. Then, the 
activation level u at the shell m should approximately satisfy the condition 
g(u) = f(u,h) + D(\ — ku) = 0. Stationary localized patterns are found 
within the parameter region where function g(u) has three real roots. The 
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boundaries of the pinning region in the plane (k, h) are determined by the 
equations (see |30|) 

u 4 + 2u (D — u 2 ) + (u 2 - D) 
~ ~Du 2 ' 

h = ^-l)+2«. (4) 

Thus, the fronts of type / are pinned within the region which lies between two 
red curves in the bifurcation diagram in fig. [2] and comprises parts b and c. 
Stationary patterns around the tree root, corresponding to the pinned front 
of type I are shown in fig. [T](a). 

Repeating the analysis for a front of type II, one can derive the pinning 
condition given in the parametric form by 

u 4 + (D + u 2 ) - 2u (D + u 2 ) 

K — 



D(l-u) 2 

h = -— ^— 2 +2u. (5) 
(1- u) 

The fronts of type II are pinned between two blue curves in fig. |2j i. e. in the 
areas c and d; an example of a stationary pattern in this case is shown in 
fig- Bb). 

Thus, our approximate theory has allowed us to identify regions in the pa- 
rameter plane (h, k) where stationary patterns are localized around the root 
or at the periphery of the trees, without the feedback control. Furthermore, 
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Figure 2: (Color online) Bifurcation diagram in the plane (h, k) for D = 0.02. 
Red curves correspond to eqs. Q, while blue curves correspond to eqs. 
The black dashed line determines the boundary on which the velocities of 
both types of fronts moving towards the tree root become equal. 
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Figure 3: (Color online) Evolution plots in absence (a-g) and in presence 
(a'-g') of global feedback are shown for different regions of the bifurcation 
diagram (fig. [2]). Global negative feedback had the intensity fx = 0.001 and 
was characterized by Sq = 100. The diffusion constant was D = 0.02. The 
degrees were (a-e) k = 4 and (f,g) k = 13. The thresholds were chosen as 
(a) h = 0.1, (b) h = 0.25, (c) h = 0.47, (d) h = 0.7, (e) h = 0.9, (f) h = 0.4 
and (g) h = 0.7. The same color coding as in fig. [l] is used. 

using fig. [2] we can consider evolution of various perturbations in different 
regions of the bifurcation diagram. 

The evolution plots in figs, pi a)— (g) show the development of the initial 
perturbation, localized in the shells between l\ = 150 and 1% = 250, in 
different regions of the bifurcation diagram in fig. [2j In absence of feedback 
control, for parameters within the region a, the activation spreads in both 
directions, towards the root and the periphery, eventually transferring the 
entire tree into the active state (fig. ^ a)). For the parameters within the 
regions b and c of the bifurcation diagram stationary patterns are formed. 
In region b, the activation spreads only towards the root; once the root gets 
activated the pattern becomes stationary and includes all nodes in the shells 
l 2 > I > (fig. |3|b)). In region c, the stationary pattern remains localized 
within the initial perturbation interval 1% > I > l\ (fig. [3](c)) . In region 
d, the activation gets annihilated, because it is pinned on the root's side 
while retreating from the periphery (fig. [3^d)). In region e, the perturbation 
retreats on both sides and disappears (fig. [3|e)). In region /, the active 
domain is gradually broadened while traveling in the root direction, but it is 
also retreated from the periphery and finally vanishes. In region g, the local 
activation is shrinking while traveling in the same direction. 
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Figure 4: (Color online) Localized stationary patterns in an ER network (a- 
d) and in a SF network (e-h) with (k) = 6 and iV = 300. The black square 
indicates the node at which the activation was initially applied. Activation 
patterns of entire networks are shown for fi = 0.004 in the upper row (a,e). 
The bottom row shows activated subnetworks at the increasing feedback 
intensities \i = 0.004 (b,f), fx = 0.008 (c,g) and [i = 0.012 (d,h). Other 
parameters are ho = 0.1 and D = 0.02; the same color coding as in fig. [I] is 
used. 



Thus, in absence of feedback control, local perturbation gives rise to sta- 
tionary patterns only for the parameters within the regions b and c of the 
bifurcation diagram in fig. [2] As we show below, in presence of global feed- 
back, the threshold h is not constant but varies with the total activation, 
according to eq. [2j Thus, different regions of the bifurcation diagram can be 
transversed, until a localized stationary pattern is established. 

When global feedback (eq. (|2])) is applied, stationary patterns develop in 
all regions of the bifurcation diagram. For the degrees k, below the inter- 
section of curves |4]) and (|5j at h = 0.5 (gray shaded region in fig. [2]), such 
patterns stay localized around the initial perturbation. When the parameter 
ho is chosen within region a, the activation starts to spread in both directions. 
As a result, the threshold h increases, thus transferring the dynamics into 
the region b, where the outer front gets pinned and spreading continues only 
inwards to the root. As h increases further, the region c is entered, where 
the pattern becomes stationary (fig. [3fa'))- In fig. [^(b'), a similar scenario 
takes place, with the feedback gradually shifting the dynamics from region b 
to c. In region c, the initial activation remains frozen (figs. J^c')). In region 
it d, the activation does not spread in the direction of the root, but it is re- 
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treated from the periphery. Therefore, h decreases and a stationary pattern 
is formed, once h enters region c (fig. |3jd')). Choosing h within region e, we 
find that (fig. j^e')) the activation first shrinks on both sides. Then, it gets 
pinned on the root's side and after that the same evolution as in fig. [^d') 
takes place. 

In the trees with the larger branching factor, i.e. the larger degrees k, 
feedback induces localized patterns which travel towards the root. Once the 
root becomes activated, the size of the patterns changes in such a way that 
the threshold h enters region b of the bifurcation diagram and the patterns 
become stationary, as shown in figs. (3^f) and (g ; ). This behavior, leading to 
pattern formation around the tree root, is also observed for the larger degrees 
k, where pinning of fronts cannot occur. 

Summarizing, we have shown that global negative feedback can induce 
stationary localized patterns in the trees, even when they do not exist in its 
absence. The pattern formation mechanism could be understood by using 
the bifurcation diagram in fig. [2] 

4 Random networks 

Effects of global negative feedback have been numerically studied for random 
ER and SF networks (fig. [4]). As we have seen, the activation that was 
applied to one node (marked by a square in fig. [4j started to spread out. 
This produced a growing cluster of activated nodes. Its growth was however 
accompanied by an increase of the negative feedback. Consequently, the 
threshold h increased with the size of the cluster, making the activation more 
difficult. As a result, the growth of the activated cluster was slowed down and 
finally stopped. Thus, a stationary activation pattern, representing a small 
subnetwork embedded in the entire network (see fig. [4]) became formed. 

By retaining only the nodes with sufficiently high activation level u > 0.7, 
active subnetworks can be identified. A sequence of such subnetworks, ob- 
tained under an increase of the global feedback intensity, is shown in figs. |4](b— 
d) for an ER network and in figs. |4]T-h) for a SF network. As the feedback 
gets stronger, active subnetworks decrease in size and approach a tree struc- 
ture. Our statistical analysis has revealed that, in both ER and SF networks, 
the average size of active subnetworks approximately followed the power law 

5 oc /U -1 . 
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Figure 5: (Color online) Stationary patterns localized on the hubs (a,c) and 
on the periphery nodes (b,d) of the ER (a,b) and the SF (c,d) networks with 
(k) = 6 and N = 300 for fi = 0.01. Dependences of the mean activation 
< u > on the degree k of the nodes in the localized patterns for the feedback 
applied only to the periphery nodes with k < k where (a) k = 9 and (c) 
ko = 6, or for the feedback applied to the hubs with the degrees exceeding 
(b) ko = 6 and (d) ko = 7. Other parameters are ho = 0.1 and D = 0.02; the 
same color coding as in fig. [I] is used. 



4.1 Local feedback schemes 

So far, we have assumed that the negative feedback was global. However, it 
is also possible that the feedback acts only on a subset of nodes. Suppose for 
example that the parameters hi, which control the dynamics of individual 
nodes i, have different sensitivity to the global control signal depending on 
the degrees of such nodes, i.e. that 

hi = ho + vB{h){S-So), (6) 

where B(ki) is the step function, B(k) = 1 if k < ko and B(k) = otherwise. 
In this case, the negative feedback is applied only to the periphery nodes 
with small degrees k, where k < ko- By introducing such local feedback, we 
can control directly the dynamics of periphery nodes, whereas hubs remain 
non-affected and have a constant threshold h = h sufficient for activation 
spreading. The simulations show that, in this situation, the activation ap- 
plied to a hub node starts to spread forming a cluster of activated nodes 
with high degrees k > ko, whose further growth to the periphery nodes is 
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suppressed by the feedback. Thus, an active stationary cluster, localized on 
the hubs and shown for an ER network in fig. [SFa), becomes formed. 

Alternatively, we can choose function B(ki) as B(k) = 1 for k > k and 
B{k) =0 otherwise, so that the feedback is applied only to the hubs. Then, 
the evolution leads to the formation of a stationary pattern, localized on the 
periphery nodes, as shown for the ER network in fig. [5^b). 

Similar behavior was observed when local feedbacks were applied to SF 
networks. Figure [5]^c) shows a stationary pattern, localized at the hubs of 
a SF network, with the negative feedback directly affecting the nodes with 
degrees smaller than ko = 6. However, if the feeback is applied to the hubs 
with the degrees higher then k Q = 7, the active pattern is localized at the 
periphery nodes, as shown in fig. |5](d). 



5 Discussion 



Our study has shown that control and purposeful design of nonequilibrium 
patterns is possible on the networks of bistable elements. Stationary localized 
patterns were found in the trees as well as in the random ER and SF net- 
works, by applying different schemes of negative feedback. They consisted of 
nodes with high activator density, surrounded by nodes where the density of 
activator was much lower, and were representing the network analogs of sta- 
tionary spots in classical react ion- diffusion media. Their size and structure 
could be controlled by varying the feedback intensity. 

In the special case of regular trees, an analytical theory could be con- 
structed, providing the pinning conditions as a function of control parameter. 
As revealed through numerical simulations, the global negative feedback can 
drive the system dynamics into the pinning regions, thus giving rise to the 
formation of stationary patterns. 

For random ER and SF networks, feedback-induced stationary patterns 
are localized on subnetworks of the entire system. The structure and the 
size of such subnetworks can be controlled by the feedback intensity. Their 
size decreases as the feedback becomes stronger and follows a power law. For 
sufficiently high feedback intensities /i, the activated subnetwork approaches 
a tree. 

For large random networks, we have also studied the effects of negative 
feedbacks which were local, i. e. acting on the subsets of nodes with small or 
large degrees. These feedback schemes resulted in stationary patterns which 
were not small subnetworks and consisted of large groups of hub or periphery 
nodes. 

Thus, the effects of different feedback schemes on dynamics of networks 



10 



formed by diffusively coupled bistable elements were investigated. In the 
future, it would be interesting to extend this analysis and consider time- 
delayed feedback schemes, which may lead to oscillating patterns or other 
complex dynamical structures on the networks. 
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